Preliminaries

We will load the necessary libraries

library(ggplot2)
library(dplyr)

## 
## Attaching package: 'dplyr'

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

library(gridExtra)

## 
## Attaching package: 'gridExtra'

## The following object is masked from 'package:dplyr':
## 
##     combine

library(plotly)

## 
## Attaching package: 'plotly'

## The following object is masked from 'package:ggplot2':
## 
##     last_plot

## The following object is masked from 'package:stats':
## 
##     filter

## The following object is masked from 'package:graphics':
## 
##     layout

library(tidyr)

Setting wd and loading databases

setwd('..')
path <- getwd()
setwd(paste(path, "/Results/PCA", sep = ""))
#Read data
coeffs   <- read.csv("coeffs.csv")
eigenvec <- read.csv("eigenvectors.csv", header=F)
eigenval <- read.csv("eigenvalues.csv", header=F)
means    <- read.csv("means.csv", header=F)
facets   <- read.csv("facets.csv", header=F)

setwd('..')
path <- getwd()
setwd(paste(path, "/Results/FSD", sep = ""))
load("fitmodel.RData")
vectors <- read.csv("vectors.csv", row.names = 1)

Data cleaning

For data cleaning, we will remove every individual that is neither Male or Female, and individuals with NAs in any variable (Height, Weight, and Age)

#Remove no sex participants
coeffs = filter(coeffs, Sex == 'Female' | Sex == 'Male')
coeffs = droplevels(coeffs)
#Removed NAs from database
coeffs = na.omit(coeffs)
dim(coeffs)

## [1] 5374   92

After cleaning, we end up with 5374 individuals

Face plots

Facial sexual dimorphism

fsdcoords <- data.frame(Tot.male = rep(0, nrow(means)), Tot.female = rep(0, nrow(means)), 
                        Allo.male = rep(0, nrow(means)), Allo.female = rep(0, nrow(means)),
                        Nallo.male = rep(0, nrow(means)), Nallo.female = rep(0, nrow(means)))
fsdcoords[,c(1,3,5)] <- apply(as.matrix(vectors[3:5,] * -2) %*% t(as.matrix(eigenvec)), 1, function(x) x + t(means))
fsdcoords[,c(2,4,6)] <- apply(as.matrix(vectors[3:5,] * 2) %*% t(as.matrix(eigenvec)), 1, function(x) x + t(means))
averagecoords <- t(( colMeans(vectors[1:2,]) %*% t(as.matrix(eigenvec))) + t(means))

source("PlotFaces.R")
Plot2Faces(fsdcoords[ ,1], fsdcoords[,2], facets, "Total SD")

Plot2Faces(fsdcoords[ ,3], fsdcoords[,4], facets, "Allometric SD")

Plot2Faces(fsdcoords[ ,5], fsdcoords[,6], facets, "Non allometric SD")

source("Distances.R")
eucdist <- data.frame(Totdist = rep(0, nrow(means)), 
                      Allodist =  rep(0, nrow(means)), 
                      Nallodist = rep(0, nrow(means)))

eucdist[,1] <- getDistance(fsdcoords[ ,1], fsdcoords[,2])
eucdist[,2] <- getDistance(fsdcoords[ ,3], fsdcoords[,4])
eucdist[,3] <- getDistance(fsdcoords[ ,5], fsdcoords[,6])

#On the same scale
eucdist <- gather(eucdist, type, dist)
eucdist$dist <- scales::rescale(eucdist$dist)

source("PlotFaces.R")
Plot1Face(averagecoords, facets, colormap = c(as.matrix(eucdist %>% filter(type == "Totdist") %>% dplyr::select(dist))), 
          title = "Total SD")

Plot1Face(averagecoords, facets, colormap = c(as.matrix(eucdist %>% filter(type == "Allodist") %>% dplyr::select(dist))), 
          title = "Allometric SD")

Plot1Face(averagecoords, facets, colormap = c(as.matrix(eucdist %>% filter(type == "Nallodist") %>% dplyr::select(dist))), 
          title = "Non allometric SD")